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^ ■ Abstract 
O " 

o 

A new method to perform numerical simulations of light-front Hamiltonians 

cn 

0^ ■ formulated on transverse lattices is introduced. The method is based on a 

DLCQ formulation for the (continuous) longitudinal directions. The hopping 
^ . term in the transverse direction introduces couplings between fields defined on 

neighboring 1 + 1-dimensional sheets. Within each sheet, the light-cone imag- 
inary time evolution operator is calculated numerically with high precision 
using DLCQ. The coupling between neighboring sheets is taken into account 
using an initial value random walk algorithm based on the ensemble projector 
Monte Carlo technique and a checkerboard decomposition for the time evolu- 
tion operator. The structure functions of theory in 2 -|- 1 dimensions are 
studied as a trial application. The calculations are performed with up to 64 
transverse lattice sites. No Tamm-Dancoff truncations are necessary. 
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I. INTRODUCTION 



Euclidean lattice QCD allows one to calculate ground state properties of liadrons, but 
real time response functions are somewhere between difficult and impossible. On the other 
hand, deep inelastic scattering gives us information about correlation functions on (or ex- 
ceedingly close to) the light-cone. Light-front quantization seems to be a promising tool to 
describe the immense wealth of experimental information about structure functions for a 
variety of reasons: (1) correlation functions along the light-cone become "static" (i.e. equal 

= (x^+x^) I -\/2) observables in this approach. (2) structure functions are easy to evaluate 
from the light-front wavefunctions. (3) these structure functions are easily interpreted as 
light-front momentum densities. However, before one can apply the light-front formalism to 
QCD and other field theories, one has to remove the divergencies first (i.e. regularize and 
renormalize). One interesting idea in this direction is called the "transverse lattice" 
Instead of discretizing all four space-time directions (like in Euclidean lattice QCD) or the 
three space directions (like in the Hamiltonian formulation of lattice QCD) one discretizes 
only the two transverse (x^ and x^) directions while leaving the longitudinal directions (x° 
and x^) continuous. On the one hand the transverse lattice thus provides a gauge invariant 
UV regularization scheme and on the other hand it is still possible to perform canonical light- 
front quantization — making it a promising approach towards performing non-perturbative 
calculations of deep inelastic structure functions. 

There remains the question what one should do with the still continuous (i.e. infinitely 
many degrees of freedom) longitudinal directions. First one may be tempted to discretize 
the x~ = (x*^ — x'^ )/V2 (the light -cone space-) direction as well. However, such an approach 
is suffering from a fundamental difficulty: the longitudinal momentum is not conserved on 
a longitudinal lattice. Due to Bragg reflections it is only conserved modulo 2p^^^ = 2n/aL, 
where is the longitudinal lattice spacing. Normally (i.e. in normal coordinates) this is 
not a problem because the minimum of the kinetic energy occurs at P = 0. However, the 
light-cone energy decreases with increasing momentum {P~ = M^/2P+, in the continuum) 



2 



i.e. a minimum is reached for F"*" = oo! Of course on a lattice the momentum cannot 
become infinite but still the minimum of the kinetic term occurs around (depending on the 
precise form of the lattice action) p+ = vtaaxl'^ — /'^cll- Since the total momentum is not 
conserved this implies that particles tend to accumulate near this minimum. However, with 
such a large momentum the particles can resolve the granular structure of the lattice and no 
meaningful continuum limit will be obtained. Note that a similar pathology would occur in 
an unconstrained (total momentum allowed to vary) variational calculation of the light-front 
energy of a hadron in the continuum. It is conceivable that adding a Lagrange multiplier 
proportional to the total light-front momentum to the lattice action cures the problem (in 
the continuum this amounts to minimizing P = + XP^ instead of P^)- However, this 
idea will not be pursued here any further. 

Instead, I found it more useful to work in momentum space as far as the longitudinal 
direction is concerned because this allows one easily to maintain longitudinal momentum 
conservation — a crucial necessity for light-front calculations as we have seen above. One 
momentum space technique which has been widely applied to light-front quantized 1 + 1- 
dimensional field theories is discrete light-cone quantization (DLCQ) [Q. There one puts 
the system into a longitudinal box and imposes periodic or antiperiodic boundary condi- 
tions (f){x~ + L) = ±0(x~). The momenta thus become discrete and solving the equations 
of motion for a fixed value of P~^ has been reduced to diagonalizing a finite matrix (note 
that all light-cone momenta are positive and thus there is only a finite number of states 
associated with a given value of The longitudinal continuum limit is achieved by mak- 
ing P+ sufficiently large. In 1 + 1 dimensions this technique was very useful and effective 
[^J^ For example, it has been used to demonstrate the existence of a nucleon-nucleon bound 
state in QCDi^i with S'f/(2)-color and S'f/(2)-flavor in a calculation on the level of quarks 
— despite the small binding energy (^ 1%) of this "deuteron" 0. One mayor obstacle for 
applications of DLCQ to 2 + 1 and 3 + 1 dimensional field theories has been the exponential 
growth of the number of basis states with the number of transverse degrees of freedom. For 
example, a DLCQ calculation for a scalar field (with antiperiodic boundary conditions in 
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the x~ direction) on a (rather modest) 4x4 transverse lattice with a longitudinal momen- 
tum P"*" = Y has already a basis size of 779022208. A more reasonable 8x8 transverse 
lattice with the same requires a basis of 6.27 • 10^^ states! These astronomical numbers 
clearly demonstrate that any direct matrix diagonalization approach or even a Lanczos type 
algorithm is doomed to fail because one is not even able to store the wavefunction in any 
available computer. Another numerical method for studying light-front Hamiltonians is the 
light-front Tamm-Dancoff approach where one imposes severe truncations of the Fock 
space. Sometimes this method is combined with DLCQ, i.e. one formulates the light-front 
Hamiltonian in the DLCQ basis but restricts {ad hoc) the Fock space to few particle states. 
However, it is not clear to what extend this truncation modifies the dynamics and whether 
effective light-front Hamiltonians can be constructed by systematically eliminating higher 
Fock states. Because of these difficulties, all numerical studies of light-front Hamiltonians 
have been restricted to 1 + 1 dimensions and/or severe Tamm-Dancoff truncations and/or 
perturbation theory. 

In this work Monte Carlo techniques will be exploited to obtain approximate ground 
state energies and structure functions for DLCQ problems in 2 + 1 and 3 + 1 dimensions on a 
transverse lattice without any constraints on the Fock space other than those resulting form 
the discreteness of the momenta. Note that the algorithm which I will present requires that 
the interaction in the transverse direction is local (nearest neighbor interaction at most). 
That is why the theory will be formulated on a transverse coordinate space lattice. The 
algorithm would not work with a momentum space lattice. To avoid obscuring the Monte 
Carlo algorithm with other difficulties the technique will be explained for a real scalar field 
with 0^ coupling in 2 + 1 dimensions, which is one of the most simple field theories one 
can formulate on a transverse lattice. The main reason to chose (f)2+i for illustrating this 
new method is that QCD is too complicated for a "first study case" and for demonstration 
purposes. Other theories one might think of, like QED^j^i of (p^j^i are not asymptotically 
free, i.e. there is no Bjorken scaling for deep inelastic structure functions. However, it should 
be emphasized that the technique is a priori applicable to any DLCQ Hamiltonian which is 
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formulated on a transverse lattice, provided the interactions satisfy locality in the transverse 
direction. 

The paper is organized as follows. First the DLCQ Hamiltonian for 0^ in 2 + 1 dimensions 
on a transverse lattice will be constructed. In section |T|, the infinitesimal light-front time 
evolution operator will be approximated using a checkerboard decomposition and a path 
integral in the Fock space of DLCQ will be used to project out the ground state for given 
quantum numbers. Finally, in section IV, the path integrals will be evaluated using an initial 
value random walk algorithm based on a variation of the ensemble projector Monte Carlo 
method. 



II. (/>|+i ON A TRANSVERSE LATTICE 



The Minkowsky action for the 02+i model, in the continuum, reads 

A 



cont. 



m 



4! 



(2.1) 



After discretization in the transverse direction one thus obtains 



\2n 



2a? 



(2.2) 



Upon rescaling = \fa(^n the interpretation of the transverse lattice action as the action of 
a "multiflavor" field theory in being the "flavor" index in this interpretation and being 
a canonical field defined in 1 + 1 dimensions) becomes evident 



Auatt. = J dx'^dx 



(V^n+l - V^n)' 



2a2 



(2.3) 



Light-front quantization of (p.3| ) is standard first one puts the system into an x -box 
of length L with antiperiodic boundary conditions (the associated zero-mode effects and 
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implications for the renormalization are discussed in Ref. 0). The canonical commutation 
relations 



[d-(pn{x),(pm{y)] 



1 g [a„(A;)e-^Pfe^ + 4,(A;)e*Pfe^" 



k=l 



k 



as well as the antiperiodic boundary condition are satisfied for 
where 



and the a„(/c) satisfy the usual commutation relations 



anik),al^{q) 



Finally one obtains for the light-front momentum operator 



n 

= T" H £ o^i{k)an{k) [k-]: 



and for the light-front energy 



L 



where 



P" = — E(^n + K + K,n+l) 



T _ ^ ai{k)an{k) 
" o 2^ i, _ 1 

^ k=l 2 



is the usual light-front kinetic energy on each site and 



(2.4) 



(2.5) 



(2.6) 



(2.7) 



(2.8) 



(2.9) 



(2.10) 



XSpj^P^ ^ : (4(^i) + (^n{ki)) («n(^2) + a„(A;2)) {ai{k^) + a„(A;3)) (a1,(A;4) + a„(A;4 



47ra4! 



fcl,fc2,A;3,fc4=l 



^4 9 



(2.11) 
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is, up to the factor 1/a the self interaction for ^p\^i. ^Pf,p, is a momentum conserving 
Kronecker 6. Neighboring sites are coupled through the hopping term 



1 



oo 



K,n+i - ^ 1. • (2-12) 



Note that, as in all DLCQ problems, the length of the box factorizes completely. Thus we 
will in the following work with the rescaled operators 

K = P+L/2'K (2.13) 

and 

H = p-2Ti/L. (2.14) 

At least in principle one could now proceed as follows: for fixed K {K and H commute) one 
diagonalizes H, yielding the eigenvalue Ei and thus the invariant masses of the physical states 
Mf = 2KEi. Physical observables like structure functions are obtained by calculating the 
appropriate matrix element in these states. The continuum limit is reached by extrapolating 
to K ^ oo. For one or two sites this can be easily done. For four sites this is also 
still possible. However, beyond about 8 sites direct diagonalization methods (including the 
Lanczos algorithm) soon become useless due to the exponential growth of the required basis 
size. What is needed at this point is a Monte Carlo method that allows one to find the ground 
state (for given value of K) of the DLCQ Hamiltonian H more efficiently. Developing and 
testing such an algorithm will be the main subject in the rest of this paper. 

Before we come to the Monte Carlo algorithm I should add a few comments. First, 
although I have derived the transverse lattice DLCQ Hamiltonian only for 02+i5 similar 
expressions can be derived for 3 + 1 dimensions and/or other field theories ll|-§|. Many 
aspects are completely general. For example, for many field theories (including QCD^^i 
|1|]) the transverse "hopping" is provided by a nearest neighbor interactions. This is crucial 
for the algorithm which I will present in this paper because (i) one can easily approximate 
the infinitesimal light-front time evolution operator as will be explained in section |T| (ii) it 



is possible to "locally update" the states in a random walk algorithm which thus provides 
a computational advantage on large transverse lattices. Second, I should discuss renor- 
malization at this point. Besides the tadpoles (which are zero in light-front quantization 
PJTO[]) there is only one divergent diagram in 02+i5 namely the setting sun diagram (Fig.^. 
This diagram leads to a divergent self mass contribution while the associated wave func- 
tion renormalization is finite. Therefore, for 02+i it will be sufficient to add an appropriate 
mass counterterm to render the theory finite in the continuum limit. Here one has to be a 
little careful in DLCQ because the longitudinal momentum of the incoming line in (Fig.|ll) 
effectively determines which higher Fock states are allowed in the intermediate state. Hence 
the self mass, and thus the required mass counterterm, for (Fig.|l|) depend on the longitudi- 
nal momentum (of course, for large longitudinal momenta, this dependence asymptotically 
disappears). For the renormalization this leaves at least two options: one possibility is to 
calculate how the infinite part of the self energy depends on the longitudinal momentum by 
evaluating (Fig.|l|) within DLCQ. Or one allows the bare mass to depend on the longitudinal 
momentum in such a way that the physical mass of the lightest state is independent of the 
longitudinal momentum. This is possible if one uses the following "sequential" procedure: 
First one solves the DLCQ Hamiltonian for K = 1/2 (which is trivial) and adjusting the 
bare mass for K = 1/2 until one satisfies the renormalization condition. Then one repeats 
the same procedure for = 3/2 while keeping the bare mass for K = 1/2 partons fixed 
and so on. In the continuum limit (large longitudinal momenta) both methods should be 
equivalent. However, I preferred to use the second method (sequential renormalization) be- 
cause it seemed to converge faster numerically in K. Furthermore, because one must verify 
numerical convergence in K, it is anyway necessary to repeat the calculation for several 
K. In addition, the computer time spent usually grows rapidly with i.e. the numerical 
effort associated with additional calculations for small values of K in the sequential renor- 
malization is negligible. Another advantage of renormalizing such that Mphys is independent 
of K is that multiparticle thresholds in two point functions appear at the correct energy 
values relative to the single particle pole. This issue has often been neglected in DLCQ 
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calculations. Furthermore, while the first method works only for 02+i because there is only 
one divergent self-energy diagram in this theory, the sequential method is non-perturbative 
and can thus be applied to the DLCQ Hamiltonian for any field theory. 

The coupling constant for 02+i receives only a finite renormalization. Nevertheless, 
in principle one should always fix the bare coupling by imposing another renormalization 
condition. Since the spectrum of consists only of the fundamental particle and its 
scattering states, and since the mass of the fundamental particle has already been used up 
to fix the bare mass term, one has to use physical observables other than the mass spectrum 
in this case. However, since this is a special feature for 02+i ^ decided not to renormalize 
the coupling at this point. In the following, all results will be quoted with the bare coupling 
constant for which they were calculated. 



III. THE MONTE CARLO PROCEDURE 

In this section, we will use the fact that repeated application of exp{—eH), where H is 
the DLCQ Hamiltonian ( |2.14| ), on any state \ip{K) > with given light-cone momentum K 
(^TT^) gives the ground state \ipo{K) > for this particular K — of course provided li'iK) > 



is not orthogonal to the ground state \iI!q{K) >. The light-front imaginary time evolution 



operator exp(—eH) can be approximated using the Trotter formula 



£3 



2^ {[Ha, [H,, Ha]] - 2 [H,, [Ha, H,]]} + ... 
^ ^-eHa/2^-eH,^~eHa/2 ^^(^^3^ (3.1) 

where the following choice is particularly useful 

n n odd 

Hb = lY.iTn + Vn)+ E Vn,n+1. (3.2) 

n n even 

This choice is motivated because Ha and H^ individually separate into sublattices which 
contain only 2 sites each 



Ha^ -(Ti+yi+T2 + y2)+14,2 + 

^(r3 + y3 + 74 + i^4) + i^3,4 +... 

Hh^ ^ (Ts + ^2 + ^3 + Fa) + 1^2,3 + 



-(T4 + V4 + T5 + y5) + V^4,5 +••• 



= -f^2,3 + -f^4,5 + ■■■ 



(3.3) 



i.e. both Ha and Hi, are sums of commuting operators and can be easily diagonalized and 
exponentiated numerically! This task is facilitated even more by the fact that ifi_2, -f^2,3 etc. 
each satisfy longitudinal momentum conservation, i.e. -^^^2,3 etc. are all block diagonal 

in the DLCQ basis, where the blocks are labeled by the sum of momenta on sites 1 and 2, 
etc. ! 

Once one has constructed exp(— siJ) one can proceed to evaluate physical observables in 
the ground state, using 

/ \ N 

<i;^{K)\H(e-^") \UK)> 



Eo{K) = lim 



N-^'oo <V'/(X)|(e-^)^|V'i(i^) > 
for the ground state energy and 

<MK)\d\MK) > _ 
<MK)\MK) > 



(3.4) 



> 



lim 



(3.5) 



for any other observable, provided < 'ilJo{K)\'ipi{K) >^ and < ■0/(-f^)|'0o(-f^) >^ 0. The 
results thus obtained still depend on e because exp(— si?) has been approximated only up 
to C(£^). The £ — > limit can be obtained by extrapolation. 

The above products are evaluated as follows. First one inserts a complete set of states 
at each £-step, e.g. 

X < ?3JV 4jv| ■ ■ ■ e-s^"!?! ki >< ?i kili/jiiK) >, (3.6) 



< 



i/jfiK)] (e li^iiK) >= ^ < -0/(^)16 2^l4iv+i hN+i >< 13N+1 hN+i\e '"'{isn ksN 



> 
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where the \ii ki > {I = 1, 3A^ + 1) are a direct product of complete sets of states at each 
site. The vector notation is used to label these states where the n-th component of k refers 
to the longitudinal momentum on the n-th site and the n-th component of i refers to the 
internal excitation level on this n-th site with a given momentum^ For example for a lattice 
with 4 sites a randomly picked state could look like this: 

k= (0,l,f,3) 

?=(1,1,3,2) (3.7) 

which would mean that site one carries longitudinal momentum, site two 1 unit, etc. while 
the internal excitation numbers are 1 (for zero momentum there is only one state — the 
vacuum), 1, 3 and 2 respectively. Of course one can chose many other ways to enumerate 
the basis but I found the above method the easiest to implement in a computer code. 

In this work, a free Fock space basis has been selected to represent the internal excitations 
at each site because structure functions are diagonal only in a free Fock space basis. However, 
many other choices are conceivable as well. In fact, in many cases it may be more efficient 
numerically to chose some kind of interacting basis on the sites instead of a free basis. 
Particularly if one is more interested in energy eigenvalues or observables other than structure 
functions. 



Several Monte Carlo techniques are available to perform the summations in Eq.(3.6). I 



have tried the projector Monte Carlo method ||TT| as well as the ensemble projector Monte 
Carlo method |T^. In the projector Monte Carlo method, in each sweep one starts from a 
state l^i ki > which is randomly picked with probability < ii ki\ilJi{K) >. In the next step 
one picks a state \i2 k2 > with probability W{ii,ki i2,^2)- The probabilities chosen in 
this work are 



^After all there is still a 1 + 1 dimensional field theory attached to each site, i.e. for each site- 
momentum many different states are possible. 

^In principle, there is little restriction in the choice of the "probabilities". However, I found 
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u/r r ^ r \ < ii+i ki+i\U\ii ki > 

W{ii,ki^ti+i,ki+i) = — — ■ ^ (3.8) 

E- fcl < « k\U\ti ki > I 

where U = exp(-|ifj for / = 3A^ or / = 3A^ + 1 and U = exp{~eHb) for / = 3A^ + 2 
, = 0, 1, 2, .. . 2 Since U factorizes into two-site sublattices, so does the transition 
probabihty W. For example, suppose (on a lattice with 4 sites) 

ki = {kf, kf, kf, kf) 

ii= {z},zf,if,if) (3.9) 

and suppose U = exp(— |ifa). Then one first selects kf_^_'^^, if_^_-^^, kf_^^ with probability 
proportional to 

'^^l+i^kl_^_i,if^i,kf_^_^\e ^^^-^lij^kljifjkf > . (3.10) 

The actual selection can be done with a Metropolis algorithm. Then one applies the same 
procedure to sites 3 and 4. Having thus "updated" the entire lattice one proceeds to the 
next e-step where U = exp{—eHf,). This time it is sites 2 and 3 as well as sites 4 and 1 
which interact with another. And so on until one reaches < i^N+i ^37V+i I • For such a "path" 
one then computes the product 

S = Sf^sN^i X S'3jv+i,3Ar X ... X 5*2,1, (3-11) 

where 

5/,3JV+i =< ^/(K)|e-5^"|4iv+i hN+i > (3.12) 



Eq.( |3.8| ) useful numerically. 
^Note that two adjacent steps with U = exp(— |ffa) can be combined into one step with U = 
exp(-ei?a). 

^Note that momentum conservation requires kl_^_i + kf_^_'^ = kf + kf and thus restricts the possible 
choices. 
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and the "scores" at each step are the actual transition amphtudes (the matrix elements of 
U) divided by the "probabilities" 

S,,,,= <V>'y'V^'^ . (3.13) 
W{ti, h ki+i) 

These products of scores ( p.ll|) are then averaged over many sweeps. Similarly one computes 
the numerators in Eqs.( |3.4| ) ( |3.5| ). For example, if one wants to calculate a structure function, 
one averages over the product of scores times the structure function measured somewhere 
near the middle of the path and divides the result by the average of the product of scores.^ 
The disadvantage of the projector Monte Carlo method is that once a small "score" 
occurs along a path, the whole product for the path will contribute negligibly to the sum 
of paths, i.e. the path integral will be dominated by a few paths which do not contain any 
small element and the relative statistical fluctuations in the path thus increase linearly with 
the number of steps in the path ||12|,|13|. Since on the other hand one is interested in making 



the path as long as possible in order to project onto the ground state while keeping e small to 
avoid systematic errors, the projector Monte Carlo method turns out to be rather inefficient. 
A variation of this method, the ensemble projector Monte Carlo method suffers less from this 
problem. There one starts form an ensemble of states at step 1. The probabilities and scores 
are calculated for each member of the ensemble in the same way as for the projector Monte 
Carlo method. However, after each e-step, when one calculates the scores, one replicates and 
deletes members of the ensemble by the following rule: first one calculates the average score 
S (ensemble average). Then one replicates each state v in the ensemble with multiplicity 



n,, = Int 



+ 'random number G (0; 1)' 



(3.14) 



where Sjy is the score for the transition to state u. This way, any path where a very small 
score occurs will most likely be eliminated (and no further computer time will be wasted on 



^Here it becomes clear why a diagonal representation is preferable for an easy computation of 
structure functions. 
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this path) while important paths branch out and contribute with multiple weight. Note that 
the size of the population is not fixed because the states are replicated with a multiplicity 
that depends on a random number. In general, the population thus fiuctuates. Sometimes, 
the population grows or shrinks slowly. In order to keep the population stable on the 
average one can counterbalance the growing/shrinking by choosing S in Eq.( |3.14j ) a little 
larger /smaller than the average score. 

For the path integral one only has to evaluate the sign of the scores 

S = Sf^3N+i X sign {S3n+i,3N x ... x 5'2,i) . (3.15) 

The absolute values of the scores along the path are already taken into account since they 
determined the multiplicities in the random walk (|3.14]) . Observables are calculated similarly 
as in the projector Monte Carlo Method except that one now averages over both the ensemble 
and the sweeps. For the efficiency of the algorithm, it is important to start the random walk 
with a good initial guess for the ground state \ilJi{K) >. For many light-front Hamiltonians 
the exact ground state wavefunction has a sizeable overlap with the valence configuration. 
As a first try, it may thus be useful to make a valence ansatz for \ilJi{K) > and \il!f{K) >. 
This will also be the choice in this work. For more complicated theories one can also try to 
improve the valence ansatz for \ilJi{K) > by including higher Fock components perturbatively. 

IV. NUMERICAL RESULTS 

The Lagrangian for (f)2+i contains two dimensionful parameters: the bare mass mo and 
the coupling constant A (which also carries dimension of mass) and hence an arbitrary 
scale A, specifying the units in which these dimensionful parameters are measured. In the 
following, the scale will be fixed by demanding that Mp^yg = A for the lightest physical 
particle. All other dimensionful quantities (e.g. A or a^^) will be measured in these units. 
After fixing the physical mass scale, which determines the bare mass, there is still (in the 
continuum limit) one dimensionless parameter left: A/47rA. In my numerical work I have 
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considered two extreme cases: A/47rA = 1 and A/47rA = 10 corresponding to intermediate 
and strong coupling respectively. The case A/47rA = 1 turned out to be rather boring 
because the structure function of the lightest physical particle is strongly dominated by a 
"bare" excitation \ipo{K) >^ J2n(^n{K)\^ >■ Deviations from this pointlike structure can 
be well approximated by summing a chain of "setting suns" (Fig Although this result 
was reproduced in the Monte Carlo calculations, it will not be discussed here any further 
because the real strength of the Monte Carlo method lies more in the nonperturbative regime. 
A/47rA = 10 will thus be chosen throughout the rest of this paper. 

The numerical calculations were done as follows: after selecting the coupling constant 
(A/47rA = 10) and choosing a transverse lattice (characterized by the spacing aA and by 
the number of sites) and after selecting a value for the "damping parameter" e, the bare 
masses were determined by imposing the renormalization condition Mphys/A = 1. Within 
the DLCQ formalism used here this yields bare masses which show some dependence on 
the longitudinal momentum. For p = 1/2 there is no interaction and Mphys{^) = A implies 
mo(|) = A^. In the next step the bare mass ?Tio(|) ^'^^ p = 3/2 is determined by requiring 
Mphysi^) = A (keeping ml{^) fixed). Then rriK^) while keeping mK^) and "^o(|) fi^^d and 
so on, up to p = Y — the largest momentum used in this work. In this fine-tuning process, 
at every longitudinal momentum, the physical mass of the lightest particle was determined 
using the ensemble projector Monte Carlo technique with an ensemble size of 500 states and 
with 10000 e-steps. For ipi and ipf a plane wave (zero transverse momentum) of bare ground 
state "mesons" was used 

\UK) >= \MK) >= > . (4.1) 

n 

This choice (which corresponds to the valence approximation in (p^) was motivated by the 
fact that the physical ground state particle in can be interpreted as a dressed single 
particle state. 
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After 15 "thermalization steps" the energy was sampled every 5th e-step Q (t o insure sta- 
tistically independent sampling) by separately taking the ensemble average of the numerator 
and the denominator in Eq.( p.4| ). The physical mass is obtained from 

M^, = 2KE, (4.2) 

where E is the average over all energy measurements (10000/5 = 2000 in the above pro- 
cedure). Typical (statistical) errors with these parameters where of the order of 1 — 3% 
for E and thus also for the bare masses. Whenever Mp^yg, evaluated using Eq. (|4.2|) devi- 
ated significantly from A (the renormalization condition) the bare masses where adjusted 
accordingly. 

This procedure was repeated for eA = .3, .15, .075 and for lattices with NgUes = 
4, 8, 16, 32 and 64 transverse sites and for transverse spacings aA = 1, |, i. It turned 
out that the bare masses depend only very weakly on e and NgUes which made the tuning 
rather easy. The tuning was facilitated even more by the fact that the path integral can still 
be summed up numerically exactly (without Monte Carlo) for lattices with four sites (and 
p < y)- Furthermore, one can estimate the difference between 4-site lattices and larger lat- 
tices by perturbative methods before one starts the nonperturbative tuning using the Monte 
Carlo. For a = 1/4 and large K the energy measurement described above became too noisy 



due to the sign problem discussed at the end of section |T|. For these cases, I used a slightly 
different algorithm: Instead of evolving the initial state for 10000 e-steps and sampling the 
energy every fifth step, the initial state was evolved only for 10 e-steps. This procedure 
was then repeated 2000 times to obtain the same statistical sample size for the energies as 
in the first method. The advantage of the second procedure is the following. Since only a 
few (small) matrix elements of exp{—eH) are negative, it is not very likely to encounter a 
negative score in a given "path" ( ^.15| ) — unless the path is very long. With this scenario 
it is clear that procedure two (starting over and over again from the same initial state) has 



^Here each application of exp{—eHa) or exp(— ei^f,) is counted as one e-step. 
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much less of a sign problem than procedure one (continued evolution). However, typically 
it takes more e-steps to project onto the ground state from \ipi{K) > ( |4.1| ) than it takes to 
get uncorrelated energy measurements. It thus depends on the concrete example which of 
the two procedures is more efficient. 

After completing the renormalization for a given set of parameters {NgUes, a, e, A) one 
can proceed to evaluate physical observables. At this point let me introduce the "structure 
functions" for scalar fields. In analogy to definitions of parton distributions in QCD one 
can introduce 

f{x) = X r ^e^-«-^" < Hp^Mminmp^) > (4.3) 

J-00 4:71 

as the light-cone momentum density of elementary quanta in the state {ipip^) >f\- The 
normalization is such that the momentum sum rule reads Jq dx xf{x) = 1. Upon discretizing 
the structure function can be expressed as (|2.5|)(|2.6|) 

fi^p) = ^"e" < HK)\ai{p)a^{pmK) > (4.4) 

^ n=l 

where Xp = {p — ^)/K, p = 1, K + ^. This expression ( [4.4|) is the form used in this work. 
What should one expect f{x) to look like for First, for trivial kinematical reasons, 

f{x) is nonzero only for x G [0; 1]. Second, since the wavefunction renormalization in 02+i 
is finite, there is a finite probability to find the physical "meson" as a bare state. Thus, in 
the continuum limit, f{x) should contain a 5-function at x = 1 with finite coefficient (in 
the structure function plots in this work this point will always be excluded because it would 
lie outside the chosen plotframe). Besides the (5-function one expects a continuum because 
the bare state can always split into three "partons" (via the 0^-interaction) which can split 
again and so on. While lowest order perturbation theory suggests a structure function that is 
peaked around x = 0{^), higher order effects (multifragmentation) will shift the maximum 
towards smaller values of x. 



^One can in fact imagine gedanken experiments that would allow one to measure f{x) but this 
point will not be discussed here any further. 
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Structure functions were evaluated by averaging over 10000 sweeps (except for a = 0.25 
where structure functions were averaged over 20000 sweeps) where again an ensemble size 
of 500 was used. Each sweep consisted of 30 e-steps. Also ipi = ipf = Eg. ([4.1[) was used 
again. The operator to measure the structure functions was inserted after the first 15 e-steps 
for each member of the ensemble.^ For the remaining steps of the random walk, a record 
is kept of the result of this structure function measurement from which each subsequent 
ensemble member has evolved [|l3l. Finally, the result for the measurement of f{x) in this 
sweep is obtained by separately evaluating the ensemble average of the numerator and the 
denominator in Eq.( p.5|) . This result is then averaged over the sweeps. In order to investigate 
convergence with respect to the number of e-steps in the final state, this number was kept 
variable (1 — 15). Typical results for such a structure function measurement, as a function 
of the number of e-steps before taking the overlap with < iljf{K)\, are shown in Fig.|^. For 
all values of e used in this work, the plateau sets in before 10 5-steps, i.e. the ground state 
expectation value of the structure function can be extracted from these results.^ There are 
several reasons for this rapid convergence. First, even for A/47rA = 10, the true ground 
state has a large overlap with the bare state > ( |4.1| ). Furthermore, any contamination of 
l^pi > with excited states is filtered out very efficiently because higher states are suppressed 
by the square of their masses 

exp{-eH) = '^\n >< n\ exp{-e—^). (4.5) 

The lowest excited state with the same quantum numbers as {ipi > is a scattering state 
consisting of 3 ground state mesons with a threshold at Mi = SMq (i.e. Mf = 9M^ !). For 
example, after 15 steps with e = 0.3, Mq = 1 and K = y the ground state is enhanced 

^Note that the structure functions in DLCQ are defined for a discrete set of points only (8 points 
for K = Hence "measuring the structure function" means evaluating it at a few points ( [4.41 ). 

^Since tpi = tpf, conclusions about convergence in the final state can also be applied to the 
convergence in the initial state as well. 
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by a factor of e~°'^/e~^'^ ~ 11 compared to the first excited state. Since most of the 
"contamination" comes not from the threshold itself but from many higher excited states, 
the filtering process is even more efficient than this numerical example illustrates. 

Since the structure functions turned out to have converged already after 10 e-steps in the 
final state, Q the measurements for /(x) with 10 — 15 final steps were then averaged. Note 
that these measurements are statistically correlated. The statistical error for the average 
result was estimated by taking the statistical error for the measurement after 10 e-steps in 
the final state (the statistical error almost does not change from step 10 to 15). 

The results are shown in Figs.|^-p. First one has to make sure that the longitudinal 
momentum K was large enough. In Fig.^ results with = | and K = ^ are compared for 
some typical choices of the other parameters. One can easily imagine that the K = ^ and 
the K = Y results lie almost on the same smooth curve, indicating that = y is large 
enough. In Fig.^ it is demonstrated how the results converged as a function of the number 
of sites. Note that (with the exception of NgUes = 2 and 4, which were done numerically 
exactly) although the same number of sweeps and size of the ensemble were used for all 
values of NgUes, the statistical error bars increase only very slowly with NsUes- The same 
is true for the CPU-time required for the Monte Carlo calculation. The reason is that the 
light-cone vacuum far away from physical particles is trivial. Thus there are no statistical 
fluctuations arising from "updating the vacuum" on huge lattices (much larger than the 
transverse size of the particles). 

The e-dependence is illustrated in Fig.||. One way to understand Fig.| is to evaluate 
the double commutator which governs the 0{e^) corrections in Eq. (|3.1j ). The rather lengthy 
expression will not be given here, but one can immediately guess the basic features. First, it 
is the hopping term Vn^n+i oc which gives rise to a nonvanishing commutator in [Ha, H^], 



Note that this means that 10 e-steps would also have been sufficient in the initial state. However, 
this was not clear before the calculations were completed. 
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which "explains" the increase of the finite e corrections with decreasing a: Fig.|^ a) vs. b) (of 
course, in order to be quantitative one has to to evaluate the matrix elements which could 
show an a dependence as well from the wavefunctions). Furthermore, the double commutator 
contains terms which spread over up to 4 transverse sites. It is thus not surprising to find 
a difference between the finite e effects on 4-site and larger lattices: Fig.^ a) vs. c). Finally 
one has to take the a limit (Fig.^. This is the most difficult part because one first 
has to make sure that everything else has converged for fixed a. As discussed above, this 
requires smaller e for smaller a (hence more steps to project onto the ground state and hence 
more computer time). The required number of lattice sites also increases with decreasing 
a. Furthermore, since the coupling between sites goes like a~^, the statistical fiuctuations 
induced by this term increase as a decreases. An additional reason for the increase of the 
statistical fiuctuations with decreasing a may lie in the choice of probabilities (p.8|) . On a 
two site lattice, the coupling between the two sites is twice as strong as between adjacent 
sites on a larger lattice (one step on the lattice in either direction leads to the same result 
on a two-site lattice with periodic boundary conditions). Thus the two site amplitudes tend 
to overemphasize highly excited states (compared to the actual physical situation on a large 
lattice). Since this "driving force" away from the ground state increases with smaller a, an 
increase in the statistical fiuctuations results. Evidently, there is still room to improve the 
algorithm chosen here. However, since the purpose of this work is only a feasibility study 
and since the most useful choice for the probabilities ( |3.8|) is expected to depend on the 
theory under consideration, this point will not be elaborated on here any further. 

In general, one would extrapolate the results first to NsUes and e ^ and then 

analyze the a-dependence of the extrapolated results. Here this is not necessary because the 
above results demonstrate that aNgUes = 4 and e = 0.075 are already close enough to the 
continuum limit for a > 0.025. Thus one can directly use these results without having to 
extrapolate. 

While the structure function still changes significantly as one goes from (a, A^^ites) = (1, 4) 
to {a, N sites) = (0.5,8), there is only a slight difference (at very small x and the slight 
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"shoulder" around x ~ 0.6) between {a^NsUes) = (0.5,8) and {a,Nsites) = (0.25, 16), i.e. for 
a < 0.5 the numerical results (Fig.|^) are almost independent of a which shows that the small 
a scaling region has been reached. Notice that the wave function renormalization in (f)2+i is 
finite and thus structure functions scale to a finite limit asQ"^ ^ oo (corresponding to a ^ 0); 
i.e. there is no logarithmic evolution for 02+i (if there were logarithmic evolution then the 
probability to find the state in the valence configuration would tend to zero as oo). 
The shape of these structure functions can be understood as follows. With lowest order 
perturbation theory (fragmentation: valence state 3 partons) one obtains a structure 
function which has a maximum near x = 1/3. Once one includes the (nonperturbative) 
interactions within the three particle sector, the structure function becomes smeared out. 
The rise at very small x can only be understood from multiple fragmentation processes. The 
"shoulder" near x ~ 0.5 — 0.6 is a nonperturbative effect. It arises because the repulsive 0^ 
interaction tends to enhance components of the wavefunction with a node in the longitudinal 
direction. This also explains why the shoulder is absent for lattices with 4 or less sites. The 
reason is that parity and Bose symmetry require that a node in the longitudinal direction 
is accompanied by a node in the transverse direction. On a small lattice, components of 
the wavefunction with a node in the transverse direction have a very large kinetic energy 
and are thus suppressed. The A^sites-dependence of the results (Fig.|D indicates that the 
physical states have a transverse extension of about 1 — 2 in the above units (i.e. 4 — 8 
lattice spacings for a = 0.25) because for lattices with a larger physical volumes, there 
is no significant volume dependence of physical observables. However, for more detailed 
information one would have to measure transverse density density correlation functions. 



V. SUMMARY AND CONCLUSIONS 

I have shown that it is perfectly feasible to perform a Monte Carlo calculation for a light- 
front Hamiltonian formulated as a DLCQ problem on a transverse lattice. The transverse 
lattice was used to separate longitudinal and transverse dynamics. The dynamics within the 
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longitudinal "sheets" attached to each transverse lattice point was solved using DLCQ and 
conventional matrix diagonalization algorithms. The transverse dynamics was then included 
by means of Monte Carlo techniques. For this purpose, the light-cone imaginary time evo- 
lution operator < il)f\exp{—NeP~)\ipi > was approximated by breaking up the light-front 
Hamiltonian P~ into two terms {P~ ~ -\- Pf), each of which contains only interac- 
tions between pairs of sites. exp{—NsP~) is then evaluated by alternating application of 
infinitesimal "evolution" -operators generated by respectively. For calculating the actual 
path-integral, a variation of the ensemble projector Monte Carlo technique was used. 
In this whole program it was crucial that the longitudinal momentum P+ was conserved, 
otherwise the light-front Hamiltonian P~ has no minimum corresponding to a physically 
meaningful particle solution. This was the major reason to use DLCQ to solve the longi- 
tudinal dynamics and not, for example, a longitudinal lattice. It was furthermore crucial 
that the transverse lattice action is local, i.e. it involves only interactions between neigh- 
boring sites. On the one hand, due to the locality, it was thus possible to perform the 
abovementioned breakup of the light front Hamiltonian P~ = P~ + P^ in such a way that 
^a/h can be written as direct sums of Hamiltonians acting on 2-sites-lattices. Thus P'j^ 
can be easily diagonalized and exponentiated numerically. The locality of the transverse 
dynamics was furthermore important when updating the states at each Monte Carlo step. 
This is because locahty of the interaction allowed to formulate the updating of the states at 
each £-step as a sequence of independent local updatings (this will also be important when 
running the algorithm on parallel computers). 

The advantages of the algorithm introduced in this paper are as follows: most impor- 
tantly, structure functions are diagonal in the DLCQ basis used here and are thus easy to 
evaluate numerically. Furthermore, since the light-front momentum P"*" is manifestly con- 
served, and since the vacuum has P"*" = 0, physical particle states are always manifestly 
orthogonal to the vacuum. Thus the Monte Carlo procedure will always converge to the 
particle solution with lowest invariant mass for that particular value of P+ (and the same 
discrete quantum numbers as the initial state). In addition, since the light-front vacuum 
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is trivial, no computer time is "wasted" to solve for the vacuum surrounding a physical 
particle while one is interested in the particle only (an annoyance for very large euclidean 
lattices). Another advantage of using the light-front Hamiltonian in the Monte Carlo pro- 
cedure is that excited states are suppressed by the square of their masses: exp{—NeP^) = 
J2n \n > exp(-A^£M2/2P+) < n| (instead of exp(-/5P°) = En l'^ > exp(-/?M„) < n\ which 
one encounters in a conventional Hamiltonian formulation). Thus fewer steps are necessary 
to filter out the ground state. The transverse lattice formulation also avoids part of the 
species doubling problem for fermions because doublers occur only for the latticized trans- 
verse directions. Thus at most four species of fermions are generated if one starts from the 
naive fermion action in 3 + 1 dimensions. Hence, by means of staggering, one can easily get 
to two light flavors of fermions . 

One of the main disadvantages of the new method are the occurrence of a few nega- 
tive matrix-elements in the infinitesimal evolution operator (although its eigenvalues are of 
course still positive). Thus negative scores resulted occasionally which slightly increased the 
fluctuations of the signal. It is expected that this problem gets worse once fermions are 
introduced because of the minus sign in exchange terms (of course this then is nothing but 
the usual sign problem for fermions). Another disadvantage is the explicit breaking of the 
symmetry between longitudinal and transverse directions which makes it difficult to recover 
full Lorentz invariance in the continuum limit. Although this was no problem in 02+i' 
will be difficult to obtain Lorentz invariant physical results in general (for theories where 
the fundamental particles carry spin). 

Numerous extensions of this work are conceivable. First of all, there is no profound 
difficulty to extend the formalism from 2+1 to 3 + 1 dimensions. Of course, in 3 + 1 
dimensions P~ = + P^^ has to be decomposed in such a way that P~i^ each can be 
written as direct sums of Hamiltonians on 1-plaquette lattices. Another difference is the 
dependence of the numerical results on the lattice spacing a: since 02+i has a finite wave 
function renormalization, scaling is exact which implies that the structure functions approach 
a finite limit for a ^ 0. Of course in a renormalizable theory (like QCD) this will not be the 
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case and the structure functions will diverge as a — > — corresponding to the logarithmic 
evolution in . However, at least in principle, this is not a problem because one can always 
perform the Monte Carlo calculations with smaller and smaller spacing a until one can 
match on to the perturbative evolution. Since scaling in sets in at moderate 

values already, there is reason to expect that this is also possible in practice (i.e. numerically 
practical). 

In this work I have investigated only the projector and the ensemble projector Monte 
Carlo method because these are rather straightforward to apply to DLCQ problems. For 
more complicated field theories, it may be necessary to use more efficient techniques, like 
guided random walks [|1^] in the Monte Carlo procedure. One could also imagine combining 
the Monte Carlo technique presented in this work with renormalization group techniques. 
On the one hand this means using perturbative renormalization group arguments to facil- 
itate the determination of the effective coupling constants in the light-front Hamiltonian. 
On the other hand one could use the Monte Carlo procedure to perform nonperturbative 
studies of renormalization group flow for light-front Hamiltonians without having to resort 
to uncontrolled truncations of the Fock space. Similarly, it is conceivable that Monte Carlo 
results are helpful in determining the effective coupling constants for the Tamm-Dancoff 
approach to solving light-front field theories. Besides structure functions, one can also use 
LFEPMC to calculate valence wavefunctions which have many interesting applications to 



various exclusive hard scattering processes [Tj]. However, the main question is whether one 
can apply the light-front ensemble projector Monte Carlo to QCDj,^x. Here the main dif- 
ficulty which remains is formulating compact (to render the transverse lattice action gauge 
invariant) QCD on a transverse lattice using DLCQ or to construct another approximation 
to QCD on a transverse lattice which is suitable for DLCQ. Once one knows the DLCQ 
Hamiltonian for QCD, it is straightforward to apply the Monte Carlo technique presented 
in this work. 
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FIGURES 

FIG. 1. Divergent self-energy diagram in ^2+1 

FIG. 2. Result of structure function measurements for 16 transverse sites, a transverse spacing 
of a = 0.5 and a longitudinal momentum oi K = 15/2 as a function of the number of e-steps before 
taking the overlap with i/jf. Because of the smaller value of £ in b), more £-steps are necessary to 
project on the ground state. 

FIG. 3. Typical examples for structure function measurements to illustrate the (in)dependence 
of the results K = 9/2 and K = 15/2. 

FIG. 4. Dependence of the structure function on the number of transverse sites for the smallest 
lattice spacing (a = 0.25) and the smallest e-step used in this work. The longitudinal momentum 
is K = 15/2. Both plots correspond to the same structure function but at different x-values. In 
order to avoid overlapping symbols from different x-values, the results are displayed in two plots. 

FIG. 5. e-dependence of the structure functions for K = 15/2. a) for 4 sites and a = 0.5; 

b) same as a) but for a = 0.25; c) same as a) but for 8 sites. Note the statistically significant 
deviation of the results for e = 0.3 in c) in the intermediate x region, while a) (same a as c) shows 
no such deviations. 

FIG. 6. a-dependence of the structure functions for e = 0.075 and K = 15/2. Note that 
aN sites — the physical volume — is kept fixed in a) , b) and c) . 
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